--- title: Data Transforms for Hyperspectral Datacubes keywords: fastai sidebar: home_sidebar summary: "Hyperspectral datacubes contain three axes (cross-track, along-track, and wavelength) and are stored in 3D numpy ndarrays. When using a pushbroom scanner, the datacube are filled gradually. The implementation here is based on a circular buffer and we provide additional methods to apply a pipeline of transforms which, for example, can be used for smile correction, and radiance conversion. " description: "Hyperspectral datacubes contain three axes (cross-track, along-track, and wavelength) and are stored in 3D numpy ndarrays. When using a pushbroom scanner, the datacube are filled gradually. The implementation here is based on a circular buffer and we provide additional methods to apply a pipeline of transforms which, for example, can be used for smile correction, and radiance conversion. " nb_path: "00_data.ipynb" ---
For example, we can write to a 1D array
cib = CircArrayBuffer(size=(7,),axis=0)
for i in range(9):
cib.put(i)
cib.show()
for i in range(9):
print(i,cib.get())
Or a 2D array
plots_list = []
cib = CircArrayBuffer(size=(4,4),axis=0)
cib.put(1) # scalars are broadcasted to a 1D array
for i in range(5):
cib.put(cib.get()+1)
plots_list.append( cib.show().opts(colorbar=True,title=f"i={i}") )
hv.Layout(plots_list).cols(3)
The OpenHSI camera has a settings dictionary which contains these fields:
camera_id is your camera name,row_slice indicates which rows are illuminated and we crop out the rest,resolution is the full pixel resolution given by the camera without cropping, andfwhm_nm specifies the size of spectral bands in nanometers,exposure_ms is the camera exposure time last used,luminance is the reference luminance to convert digital numbers to luminance,longitude is the longitude degrees east,latitude is the latitude degrees north,datetime_str is the UTC time at time of data collection,altitude is the altitude above sea level (assuming target is at sea level) measured in km,radiosonde_station_num is the station number from http://weather.uwyo.edu/upperair/sounding.html,radiosonde_region is the region code from http://weather.uwyo.edu/upperair/sounding.html, andsixs_path is the path to the 6SV executable.The pickle file is a dictionary with these fields:
camera_id is your camera name,HgAr_pic is a picture of a mercury argon lamp's spectral lines for wavelength calibration,flat_field_pic is a picture of a well lit for calculating the illuminated area,smile_shifts is an array of pixel shifts needed to correct for smile error,wavelengths_linear is an array of wavelengths after linear interpolation,wavelengths is an array of wavelengths after cubic interpolation,rad_ref is a 4D datacube with coordinates of cross-track, wavelength, exposure, and luminance,sfit is the spline fit function from the integrating sphere calibration, andrad_fit is the interpolated function of the expected radiance at sensor computed using 6SV.These files are unique to each OpenHSI camera.
For example, the contents of CameraProperties consists of two dictionaries. To produce the files cam_settings.json and cam_calibration.pkl, follow the steps outlined in the calibration module.
cam_prop = CameraProperties("assets/cam_settings.json","assets/cam_calibration.pkl")
cam_prop
cam_prop.calibration["rad_ref"]
We can apply a number of transforms to the camera's raw data and these tranforms are used to modify the processing level during data collection. For example, we can perform a fast smile correction and wavelength binning during operation. With more processing, this is easily extended to obtain radiance and reflectance.
Some transforms require some setup which is done using CameraProperties.tfm_setup. This method also allows one to tack on an additional setup function with the argument more_setup which takes in any callable which can mutate the CameraProperties class.
You can add your own tranform by monkey patching the CameraProperties class.
@patch
def identity(self:CameraProperties, x:np.ndarray) -> np.ndarray:
"""The identity tranform"""
return x
If you don't require any camera settings or calibration files, a valid transform can be any Callable that takens in a 2D np.ndarray and returns a 2D np.ndarray.
Depending on the level of processing that one wants to do real-time, a number of transforms need to be composed in sequential order. To make this easy to customise, you can use the pipeline method and pass in a raw camera frame and an ordered list of transforms.
To make the transforms pipeline easy to use and customise, you can use the CameraProperties.set_processing_lvl method.
cam_prop = CameraProperties("assets/cam_settings.json","assets/cam_calibration.pkl")
cam_prop.set_processing_lvl(-1) # raw digital numbers
test_eq( (772,2064), np.shape( cam_prop.pipeline(cam_prop.calibration["HgAr_pic"])) )
cam_prop.set_processing_lvl(0) # cropped
test_eq( (452,2064), np.shape( cam_prop.pipeline(cam_prop.calibration["HgAr_pic"])) )
cam_prop.set_processing_lvl(1) # fast smile corrected
test_eq( (452,2054), np.shape( cam_prop.pipeline(cam_prop.calibration["HgAr_pic"])) )
cam_prop.set_processing_lvl(2) # fast binned
test_eq( (452,108), np.shape( cam_prop.pipeline(cam_prop.calibration["HgAr_pic"])) )
cam_prop.set_processing_lvl(4) # radiance
test_eq( (452,108), np.shape( cam_prop.pipeline(cam_prop.calibration["HgAr_pic"])) )
cam_prop.set_processing_lvl(6) # reflectance
test_eq( (452,108), np.shape( cam_prop.pipeline(cam_prop.calibration["HgAr_pic"])) )
cam_prop.set_processing_lvl(5) # radiance conversion moved earlier in pipeline
test_eq( (452,108), np.shape( cam_prop.pipeline(cam_prop.calibration["HgAr_pic"])) )
cam_prop.set_processing_lvl(3) # slow binned
test_eq( (452,108), np.shape( cam_prop.pipeline(cam_prop.calibration["HgAr_pic"])) )
DataCube takes a line with coordinates of wavelength (x-axis) against cross-track (y-axis), and stores the smile corrected version in its CircArrayBuffer.
For facilitate near real-timing processing, a fast smile correction procedure is used. An option to use a fast binning procedure is also available. When using these two procedures, the overhead is roughly 2 ms on a Jetson board.
Instead of preallocating another buffer for another data collect, one can use the circular nature of the DataCube and use the internal buffer again without modification - just use DataCube.push like normal.
All data buffers are preallocated so it's no secret that hyperspectral datacubes are memory hungry. For reference:
| along-track pixels | wavelength binning | RAM needed | time to collect at 10 ms exposure | time to save to SSD |
|---|---|---|---|---|
| 4096 | 4 nm | ≈ 800 MB | ≈ 55 s | ≈ 3 s |
| 1024 | no binning | ≈ 4 GB | ≈ 14 s | ≈ 15 s |
In reality, it is very difficult to work with raw data without binning due to massive RAM usage and extended times to save the NetCDF file to disk which hinders making real-time analysis. The frame rate (at 10 s exposure) with binning drops the frame rate to from 90 fps to 75 fps. In our experimentation, using a SSD mounted into a M.2 slot on a Jetson board provided the fastest experience. When using other development boards such as a Raspberry Pi 4, the USB 3.0 port is recommended over the USB 2.0 port.
timestamps = DateTimeBuffer(n=8)
for i in range(8):
timestamps.update()
timestamps.data
For example, if this data was collected on 2021/10/03 at 12:18:30.010568 UTC, then there will be two files saved under the directory 2021_10_03. The two files are the NetCDF file 2021_10_03-12_18_30.nc and PNG image 2021_10_03-12_18_30.png. You can also had a custom prefix and suffix to the file name.
The plot_lib argument hs
Choose matplotlib if you want to use
Choose bokeh if you want to compose plots together and use interactive tools.
n = 256
dc = DataCube(n_lines=n,processing_lvl=2,json_path="assets/cam_settings.json",pkl_path="assets/cam_calibration.pkl")
for i in range(200):
dc.put( np.random.randint(0,255,dc.settings["resolution"]) )
dc.show("bokeh")
dc.save("assets")
import shutil, os
shutil.rmtree("assets/"+[ f for f in os.listdir("assets/") if "20" in f and "_" in f ][1])
dc.nc
import os
folder = "assets/"+[ f for f in os.listdir("assets/") if "20" in f and "_" in f ][0]
fpath = [ folder+"/"+n for n in os.listdir(folder) if ".nc" in n][0]
dc.load_nc(fpath)
dc.show("bokeh")